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The near-field interaction between two neighboring particles is known to produce enhanced ra- 
diative heat transfer. We advance in the understanding of this phenomenon by including the full 
electromagnetic particle response, heat exchange with the environment, and important radiative 
corrections both in the distance dependence of the fields and in the particle absorption coefficients. 
We find that crossed terms of electric and magnetic interactions dominate the transfer rate between 
gold and SiC particles, whereas radiative corrections reduce it by several orders of magnitude even 
at small separations. Radiation away from the dimer can be strongly suppressed or enhanced at 
low and high temperatures, respectively. These effects must be taken into account for an accurate 
description of radiative heat transfer in nanostructured environments. 
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I. INTRODUCTION 

Blackbody radiation mediates heat exchange between bodies placed in vacuum and separated by large distances 
d compared to the thermal wavelength At = lnhcjkBT . For parallel plates, this leads to a radiative heat transfer 
(RHT) rate independent of d. However, when d <C At, the rate is enhanced by several orders of magnitude due to the 
involvement of evanescent waves. Pioneering measurements [TJ[2] revealed this phenomenon, which was first explained 
in terms of near-field fluctuations [3]. After a long series of experimental [HE] and theoretical [31I6HT4] studies, recent 
observations have accurately confirmed a l/d 2 dependence for sapphire plates at room temperature down to d ~ 1 jiva 
[5] , and also a l/d dependence for large silica spheres placed near a silica plate down to d ~ 30 nm [15] , although these 
laws can be substantially corrected by nonlocal [11], phonon [16] [17], and photonic crystal [18] effects. In this context, 
the interaction of a particle with a plate has been explored both from experimental [4j [15] [19] [20] and theoretical 
P~2] [16] [21-24 fronts. Modeling heat exchange between two [25H3T] or more [32] particles has been the subject of 
intense activity as well. 

Magnetic polarization has been claimed to dominate RHT between metallic nanoparticles [27] . However, electro- 
magnetic crossed terms (EMCTs, i.e., terms mixing the electric and magnetic particle responses) have been ignored 
so far, although they could play a leading role in transfers within heterogeneous structures. Likewise, radiative cor- 
rections in the absorption of dielectric particles deserve further consideration, as we show below. Thus, the current 
level of understanding of RHT between two particles appears to be incomplete. 

Here, we formulate a complete solution of RHT between two nanoparticles within the assumption of dipolar response. 
We show that EMCTs are dominant in combinations of metallic and dielectric particles, such as gold and SiC. We 
introduce a relevant retardation correction beyond the customary treatment of polarization fluctuations, which results 
in a sizable reduction in the predicted transfer rate. Furthermore, we show that heat losses into the environment can 
be either dominant or negligible depending on the temperature and particle composition. An accurate description of 
RHT in nanostrutured environments requires incorporating these effects, for which the two-particle system discussed 
here provides a tutorial approach as well as an estimate of the importance of EMCTs, radiative corrections, and 
interaction with the environment. 



II. DESCRIPTION OF THE MODEL 

We consider two spherical particles of radius R at temperatures T\ and T2 separated by a center-to-center distance 
d along the z direction and placed in a vacuum at temperature Tq, as shown in Fig. [I] We focus on small particles such 
that R <C <i, At, so that their responses can be described through the polarizabilities a\ and ql2 [26] (see Appendix, 
Fig. 6). RHT between the particles and the environment is produced by fluctuations in the vacuum electromagnetic 
field and the particle dipoles. We simplify the notation by combining electric and magnetic field components acting 



FIG. 1: (color online). Description of the system under study. Two particles of radius R at temperatures T\ and T2 separated 
by a distance d are placed in a vacuum at temperature To . Each particle exchanges thermal energy with the other particle and 
with the surrounding vacuum. 
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on each particle j = 1, 2, as well as electric (p) and magnetic (m) dipoles, in the vectors 
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respectively. Likewise, the polarizability tensor becomes 



afh 




afh 



where I3 is the 3x3 identity matrix and the E (M) subscript refers to electric (magnetic) components. 
The net power absorbed by particle 1 is the sum of dipole and field fluctuation terms, 



p 1 = pfi eld _|_ pdi P> 



More precisely (see Appendix), 



field 
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(1) 



represents the work exerted by the fluctuating field on particle 1. Here, k = uj/c is the wave vector of light at frequency 
uj. Moreover, the self-consistent fields Ej include the response of the system to the fluctuating source fields E^ via 
the relations 

E x =Ef + g 12 a 2 E 2 , 
E 2 = El + g 21 a 1 E 1 , 

where Qi 2 is the distance-dependent dipole-dipole inter-particle interaction, 

(A -C \ 
A C 
c 005000 
yi2 ~ C A 
-C A 
\0 000 05/ 



A = exp(zfcd) (k 2 /d + ik/d 2 - 1/d 3 ), 
B = exp(ikd) 2(l/d 3 - ik/d 2 ), 
C = exp(ikd) (k 2 /d + ik/d 2 ), 

and Q 2 \ takes the same form as Q\ 2 with C replaced by —C. In Eq. ([!]), () represents the average over field fluctuations, 
which we perform by applying the fluctuation-dissipation theorem (FDT) (33[|34] (see Appendix) 



n (uj) 



where rj 12 = (1/2) [Im{^i 2 + G21} + iRe{Q 12 - S21}], V21 = Vi 2 , Vn = V22 = (2/c 3 /3)I 6 , and n (uo) = [exp(fkj/k B T ) - 
is the Bose-Einstein distribution at the vacuum temperature Tq. 
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Likewise, the contribution of fluctuating dipoles is (see Appendix) 



^dip 
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-i(uj—uj )t 



(2) 



[iJ {pi (w')C/i2(w)p 2 (w)) 



2u A 
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where the first term inside the square brackets accounts for the effect of the field produced by particle 2 on particle 
1, while the second term describes the interaction of the particle dipole with the vacuum. The self-consistent dipoles 
satisfy the relations 

Pi =Pi +aiGi2P2, 

P2 =P2 +®<2G2lPl, 



where p^ is the fluctuating source dipole at particle j. The relevant FDT now becomes 



(see Appendix) 



[p>)]%* (^)) = ^hS(u-u f )S jrXj (u) 



where we use 



Xi =Im{a i }-(2fc 3 /3)|a i | 2 I 6) 



(3) 



rather than Im{aj} in order to prevent non-absorbing particles from undergoing unphysical fluctuations. We set the 
polarizabiltiy to = (3/2/c 3 )^^ (y = E,M), where t v i X is the dipole Mie scattering coefficient (see Appendix). This 
definition of aj complies with the optical theorem condition [35] Im{<x/} > (2/c 3 /3)|<x/ | 2 , where the equality applies 
to non-absorbing particles (xj = 0)- Incidentally, dipole and field fluctuations originate in different physical systems, 
and therefore, there are not crossed terms between the two of them. 

Finally, using the FDT to evaluate the integrals of Eqs. ([!]) and ([2|, we find, after some lengthy but straightforward 
algebra (see Appendix), 
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and u' = M (E) and ^ +1 (-1) when v — E (M). 

Neglecting the magnetic response (a^ = 0, C = 0), multiple scattering (^ = ^ = s = 1, ^ = 0), and radiative 
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FIG. 2: (color online). Dependence of the radiative heat transfer coefficient (HTC) on particle separation d. The HTC from 
particle 2 (right) to particle 1 (left) is defined as Vi/ST with T\ —Tq — T and T 2 =T + 5T (see Fig.fl]). (a) HTC at T = 300 K 
obtained from our full theory (solid curves), compared to the result of neglecting multiple scattering (circles), retardation effects 
everywhere (dashed curves), retardation effects except in \j (triangles), or EMCTs (dotted curves). Gold-gold (red curves) 
and gold-SiC (blue curves) dimers are considered. The particles radius is R = 500 nm. (b) Heat power lost by particle 2 (solid 
curves) and fraction of that power absorbed by particle 1 (broken curves, right scale). (c,d) Electric-electric (black), electric- 
magnetic (red), magnetic-magnetic (blue), and magnetic-electric (green) partial contributions to the HTC in the gold-gold (c) 
and gold-SiC (d) dimers. 



corrections in the particles response (xj — I m { a j})? the above expressions reduce to 

7T 

+ Re{a^B}lm{B} (n - m) 



Pi 



J udulm{af} j[/c 3 + 2Re{afA}lm{A} 



Im{af}(|A| 2 + |5| 2 /2) (ra 2 -m)}, 



where the n 2 — rt\ term describes direct RHT between the dimer particles and coincides with a previously reported 
expression [22j [25j [30] . The remaining no —n\ term accounts for heat exchange between particle 1 and the surrounded 
vacuum, partially assisted by the presence of particle 2. 
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III. RESULTS AND DISCUSSION 

We study in Figs. [2j [3ja), and [4] the heat transfer coefficient (HTC) between two particles in a dimer when particle 
1 is at the same temperature as the environment (T\ = To = T) and particle 2 is at a slightly different temperature 
(T2 = T + 5T). The HTC to particle 1 is defined per unit of temperature difference as Vi/ST. Under these conditions, 
only the terms Y v 2 ± and Y v 2 y contribute to the transfer. The results obtained from the above formalism (solid curves) 
are compared to several approximations consisting of neglecting multiple scattering between the particles (circles, 
calculated for u v = t v — s = 1, w v = 0), retardation effects everywhere (dashed curves, k = 0), retardation effects 
except in the particle response Xj defined by Eq. ([3| (triangles), or EMCTs (dotted curves). 

The distance dependence of the HTC is analyzed in Fig. [2] for a homogeneous gold dimer and for a dimer formed 
by gold and SiC particles. As a first observation, we note that multiple scattering events can be safely neglected 
in all cases. In contrast, retardation causes a dramatic boost in the HTC, which increases with particle distance 
in both types of dimers. Additionally, the \a\ 2 term of Eq. ^ contributes with a uniform decreasing factor in SiC 
particles (see Appendix, Fig. 5). (Notice that the absorption cross section is proportional to Xji whereas Im{<x/} 
describes absorption plus scattering, so the \a\ 2 term in Xj is removing scattering strength that is not associated 
with absorption.) Finally, EMCTs terms introduce additional channels of inter-particle interaction, thus resulting 
in higher transfer rates compared to the diagonal approximation (consisting of only including electric-electric and 
magnetic-magnetic terms), particularly in the heterogeneous dimer [Fig. Qa)]. Figure [2^c) clearly shows that the 
magnetic-magnetic terms are dominant in the homogeneous gold dimer, in agreement with previous predictions [27] , 
because metallic particles mainly contribute through magnetic polarization. This is unlike the heterogeneous cluster, 
in which magnetic-electric terms are dominant [Figure J2^d)] , thus picking up a dominant electric polarization from 
the SiC particle. Incidentally, the HTCs from gold to SiC and from SiC to gold are nearly identical (see Appendix, 
Fig. 8). 

Similar conclusions are extracted from the temperature dependence of the HTC, represented in Fig.|3ja) for a small 
distance d = 2 /im<C At = 14 /im— 14mm. Notice however the dramatic reduction in the transfer rate produced by 
retardation at high temperatures. 

The full dependence on the particle temperatures for a vacuum at To = is studied in Fig. j3^c,d). Interestingly, 
particle 1 gets cooled down (blue regions) even if particle 2 is at a higher temperature. This is due to radiation losses 
into the vacuum. However, particle 2 in the homogeneous cluster is rather efficient in transferring energy to particle 
1 and compensating for radiation losses, so that the curve separating gains (red) from losses (blue) is closer to the 
Ti = T 2 line (dashed) in that dimer [Fig. ^c)]. 

It is useful to analyze the spectral contribution of different photon energies to the HTC. At a low temperature 
T = 10 K [Fig. |4^a)], the exchange is dominated by low photon energies, for which the particles polarization show a 
featureless behavior and Xj almost coincides with Im{aj} for the size of the particles under discussion (see Appendix, 
Fig. 5). At high temperature T = 1000 K [Fig. Qb)], optical phonons emerge as a sharp infrared (IR) feature in SiC 
and plasmons show up as a broader near-IR feature in gold particles. Retardation effects also increase with T, as the 
particles appear to be large in front of At- 

An important ingredient that is often overlooked in the analysis of heat transfer relates to how much energy is 
emitted into the surrounding vacuum. We analyze this in Figs.|2|b) and[3^b) by calculating the power escaping from 
a hotter particle 2. The calculation is done by reversing the particle labels, so that only the terms j_ and n 
contribute in this case [see Eq. Q]. For the temperature of Fig. [2J just a small amount of the energy emanating 
from particle 2 ends up in particle 1. However, this fraction increases at lower temperatures [Figs. J^b)], until nearly 
complete heat transfer takes place below ~ 10 K in the homogeneous gold dimer. The fraction of heat transfer between 
the particles is thus very sensitive to temperature and particle distance (see Appendix, Fig. 9). 

IV. CONCLUDING REMARKS 

Heat dissipation in nanostructured devices is becoming a limiting factor in the design of microchips and is expected 
to play a major role in nanoelectronics, nanophotonics, and photovoltaics. Radiative losses provide a convenient 
way of handling the excess of heat produced in these devices [36] . In this context, crossed magnetic-electric terms 
and radiative corrections as those described here produce modifications in the transfer rate by up to several orders 
of magnitude, which cannot be overlooked. An analysis of how much heat is released from a dimer into a cooler 
vacuum reveals a large dependence on composition and temperature (it is strongly suppressed at low temperatures 
and dominant in hot environments). Our results for gold-SiC dimers suggest the experimental exploration of these 
effects via, for example, in-vacuum particle levitation, or by attaching one of the particles to a nanoscale tip and 
the other one to an insulating substrate. As an interesting direction, we note that RHT can be strongly modified 
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FIG. 3: (color online). Temperature dependence of RHT. (a) HTC for the particles of Fig. [5] (b) Heat power lost by the right 
particle and fraction absorbed by the left particle, (c) Power absorbed (red-yellow scale, positive) or emitted (blue-white scale, 
negative) by the left particle of a gold-gold dimer as a function of T\ and T2 with the vacuum at To = (see Fig.[T]). (d) Same 
as (c) for a gold-SiC dimer. The particle distance is d = 2 /xm and the radius is R = 500 nm in all cases. 



by the presence of additional mirrors and dielectrics that distort the exchanged electromagnetic fields. We further 
suggest the possibility of molding RHT down to the quantum regime by placing the particles in a resonant cavity. 
This directly connects to the proposed quantization of RHT [37] , similar to that observed in the conventional thermal 
conductance of narrow bridges [38] . 
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FIG. 4: (color online). Spectral dependence of the radiative HTC under the conditions of Fig.|3ja) for two different temperatures, 
as shown by labels. 



Appendix A: Particle polarizability 

We obtain the polarizability of the homogeneous spherical particles under consideration from their dipolar Mie 
scattering coefficients as = (3/2fc 3 )£!- ?1 (y = E,M), where k = lj/c. This procedure automatically incorporates a 
number of retardation corrections in the polarizability. The Mie coefficients are given by the analytical expressions 

EH 

t M = -Jl(Po)Plfl(Pl) + Pofl(Po)jl(Pl) 

1 hl +) (po)pdl(Pi) ~ Polh^MYMpiy 
t E = -3i(po)[pdi(pi)Y + e\ftoji(Po)}'ji(pi) 
1 fe, (+) (po)[piii(pi)]' - eipoh^ipo)}^^)' 

where po — kR, pi = kR^fe with Im{pi} > 0, R is the particle radius, e is its dielectric function, ji and h\ + ^ are 
spherical Bessel and Hankel functions, and the prime denotes differentiation with respect to po and p\. 

Figure [5l shows the imaginary part of the particle polarizability Im{a} compared to the absorption factor x = 
Im{a} — 2kr /3\a\ 2 for gold and SiC spheres of radius R similar to those presented in Sec. III. As expected, x approaches 
Im{a} in the kR <C 1 limit (i.e., when retardation is negligible). However, Im{a} and x behave increasingly different 
as the energy goes up. The difference between these two functions increases with particle size, or equivalently, the 
threshold for retardation effects uj = c/R (e.g., 0.4 eV for R = 500 nm) is lowered. Interestingly, the magnetic (electric) 
component is dominant in gold (SiC) particles at energies below ~ 0.1 eV, while plasmon (phonon) resonances take 
over for larger uj. 
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FIG. 5: Imaginary part of the polarizability Imla^} (solid curves) and absorption coefficient — Imja^} — (2/c 3 /3)|a^| 2 
(dashed curves) for gold (a)-(c) and SiC (d)-(f) spherical particles of different radius R, as a function of photon energy. Red 
and blue curves represent electric (y = E) and magnetic (y = M) components. The upper color scale shows the equivalent 
temperature hoj/ks- 



Appendix B: Validity of the dipolar approximation 

The dipolar approximation should be accurate in the limit of small particles compared to the light wavelength. The 
validity of this approximation is tested in Fig. [6] by comparing the extinction cross section of gold dimers calculated 
by representing the particles as dipoles (broken curves) or by including all multipoles (solid curves). The calculations 
are performed using a multiple elastic scattering of multipolar expansions (MESME) method [39]. We obtain similar 
values of the cross section in both calculations for different particle separations down to d = 1.5 fim when the photon 
energy is smaller than hc/R ~ 0.4eV, where R = 500 nm is the particle radius. The complex plasmon features showing 
up at larger energies in the full calculation are not captured within the dipolar approximation. Our formalism is thus 
appropriate to describe dimers within this range of sizes and separations for temperatures below ~ hc/kBR ~ 4580 K. 

Appendix C: Derivation of Eq. (1) 

We consider a spherical particle placed at the origin and of radius much smaller than the thermal wavelength 
At = Inkc/hsT [see Fig- 0a)], so that it can be described by its electric and magnetic polarizabilities, a E and a M , 
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FIG. 6: Extinction cross-section spectra of gold dimers calculated in the dipolar approximation (I = 1, dashed curves) and 
with full inclusion of all multipoles (solid curves) for different distances between particle centers d, as indicated by labels. The 
particles diameter is 1 fim. The cross section is normalized to the projected area of one particle. 




FIG. 7: Description of the geometry used in the derivation of Eqs. (1) and (2). 



respectively. Under illumination by an external electromagnetic plane wave of frequency uj and wave vector k, the 
total field (external plus induced) can be written 



E(r,t) =E e^ r e 



kr — iu)t 



H(r,t) =H e* k - r e 



k-r — iojt 



1 

ik 



V x £(r,cj)a M (u;)Hoe 



c.c, 



-iujt 



ttV x Q(r, oj)a E (oj)E e- ibjt + c.c, 

IK 



(Cla) 
(Clb) 



where 



g(r,u) = [fc 2 + V®V] 



[{k 2 r 2 + ikr - 1) - (k 2 r 2 + 3ikr - 3)r <g) f] , 



(C2) 



is the electromagnetic Green tensor. The work exerted on the particle by the external electromagnetic field can be 
obtained from the Poynting vector flux across a spherical surface S of radius r centered at the origin as 



47r Js 



(C3) 



Now, inserting Eqs. (CI) into Eq. (C3), we readily obtain the expression 



Pf eld = to 



2k 3 



-ta 



|Eor+w 



-la 



M 



2fc 3 



Ml 



|H | J +c.c, 



(C4) 



11 



where we have only retained stationary terms. Equation (C4) is the monochromatic version of Eq. (1). 



Appendix D: Derivation of Eq. (2) 



We now consider two neighboring spherical particles placed at the origin and at a distance d from the origin along 
the z axis, respectively, as shown in Fig.[7Jb). We assume small particles compared to both At and d, so that we can 
describe them as point electric and magnetic dipoles, and nij, respectively, where j — 1,2 labels the particles. We 
calculate the power absorbed by particle 1 from the Poynting vector flux in the direction towards the particle across 
two infinite planes A and B normal to z, as shown in Fig. [7^b). More precisely, 



-pdip _ 

1 4tt 



Ua ~ L) dK [E(r ' i) x H(r ' t)] ' (D1) 



Taking advantage of this geometry, the electromagnetic field created by the electric and magnetic particle dipoles can 
be written 

E(r, t) = !^! J2 I — [(k 2 + V <g> V) p, - iifeV x m t ] e ^ e id\z- Zi \ + c . c , (D2a) 

27F i=l,2 J Q 

H(r, t) J±-^L I — [{k 2 + V V) m, + ikV x Pi ] e ^ e ^-^\ + c .c, (D2b) 



= 1,2 

where we have used the relation 

Q ikr 



= j_ P dQ iQ.n iq \ z - Zi \ 
27rJ q 



the integrals are extended over wave vectors along x and y d irec tions, and q = \/k 2 — Q 2 (q = i\/Q 2 — k 2 ) is the wave 
vector along z for k > Q (k < Q). Finally, inserting Eqs. (D2) into Eq. (Dl) we obtain, after some straightforward 
algebra, 



V± ip = iu)A(d) (pl x P2x + Pi y P2y + rn* lx m 2x + m\ y m 2 y) + iuB(d) (p* lz P2z + rnl z m 2z ) 

2a; 4 



^dip 

— ouu 

+ iujC(d) (-p* lx m 2y +pl y m 2x +m$ x P2 V ~ m\ y p 2x ) - ^ (jp-^ 2 + |mi|" ) +c.c, 
which is the monochromatic version of Eq. (2). 



Appendix E: Fluctuation-dissipation theorem 



The fluctuation-dissipation theorem (FDT) was first formulated by Nyquist [40] and later proved by Callen and 
Welton [41]. It relates the fluctuations of the product of two operators to the dissipation expressed as the imaginary 
part of their response function. Next, we reproduce for convenience a simple derivation given elsewhere [34] with a 
notation suited to the analysis presented in this paper. 

Consider a perturbation Hamiltonian 

H'(t) = -<p{t)q(t), 

where ip(t) is a time-dependent function, q(t) is an operator in the Heisenberg picture, related to its Schrodinger 
representation qs through q(t) = ex.p(iHQt/h)qs ex.p(—iH t/h), and i^o is the unperturbed Hamiltonian. In the 
Heisenberg representation, Schrodinger's equation becomes H'\(j)) = ihd\<fr)/dt, and we have \<p) = ex.p(iHot/H)\<ps)- 
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Under the condition H'(-oo) = 0, the eigenstates of the total Hamiltonian become 



\4>m(t)) = \m) 
~ \m) 



dt'H'it')^')) 



dt'H'{t')\m), 



where the last line corresponds to first-order perturbation theory, and \m) is a state of the unperturbed Hamiltonian 
with energy E m (i.e., H \m) = E m \m)). 

The expected value of another operator p(t) is simply given by 



<P(*)> = \ Y. e ~ Em ' kBT ^{t)W)\4>m{t)) 



-E m /k B T 



(m\p(t)\m) 



dt'ip(t'){m\\p(t),q(t')}H 



(El) 



where 



-E m /k B T 



is the partition func tion at temperature T. The first term in Eq. (El) reduces to (m\p(t)\m) = (m\ps\m), and from 
here, we can recast (El) as 



(Sp(t)) = (p(t) -p(-oo)) = | dt'x(t-t>(*')> 



where 



X (t - f ) = -Mt - f )4 E e- £ - /fcBT <m| [p(t), g(f ) 



is a susceptibility. Now, using the closure relation |n)(n| = I, we can write 

(m\\p(t)A(t')]\m) = £ [<ro|p s |»> <n|g s |m> e *(^-^)(*-*')/ft _ ( m |fo|„) <n|p s |m) e -i(^-s»)(*-*')/»' 



(E2) 



Inserting this back into Eq. (|E2j), taking the time Fourier transform to work in frequency space, and using the identity 



find 



= J dt X (t) 



dt e l 



Aujt 



A + i0+ ' 



^(ra|j$s|n) (n|g 5 |ra) 



e -Em/k B T _ e -E n /k B T 

Huj + E m - E n + i0+ 



Incidentally, the zero-temperature susceptibility reads 

(0\p s \m) (m\q s \0) (0\q S \m) (m\p s \0) 



E 



_ficj + £ - + &0+ hu + E m - E + i0+ _ 
Finally, the dissipation associated with x can be written 

Im {xH} = (l - e-"»/* BT ) | E e- £ ™ /feBT <m|p s |n> <n|fo|m> <5(^ + £ m - (E3) 
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Similarly, we can write the average over fluctuations as 

S(t-f) = <p(t)«(f)> = ±J2 e ~ Em/kBT e^-^-^Hpsln) (n\q s \m) = f ^S^e 

m,n 

where 



-iuj{t — t') 



S(oj) 



2irh 



-E m /k B T 



(m\ps\n) (n\qs\m) S(hu + E m - E n ). 



(E4) 



The relation between S(uj) and Im {x(co>)} that one obtains by comparing Eqs. (E3) and (E4) constitutes the general 
form of the fluctuation-dissipation theorem: 



where 



S(u) = 2h [n(w) + l] Im{xH}, 



n(uj) 



e hou/k B T _ l ' 



(E5) 



(E6) 



is the Bose-Einstein distribution function. 

We can formulate a more useful relation by noticing that (p(t)q(t f )) is a function of t — t f , so that its double Fourier 
transform satisfies 



<pM«(o/)> = I 



= / dtdt 



I iojt+ioj't' 



s(t - 1 



') = J dr e^ T S(t) J dt' e i(u+0J '^' = 2tt 5(uj +u') S(u), 



and from here we find the expression 

(p(oj)q(oj')) = 4irh [n(u) + 1] Im{x(w)} S(oj + uj'). 



(E7) 



Proceeding as above, the Fourier transform of the fluctuation (q(t')p(t)) reads (q(oo f )p(uj)) = ex.p(—hu/kBT)S(uj), 
which together with ( E5 ) leads to 



(q{uo')p{uj)) = 47rhn(uj)lm{x(oj)} 5(uj + uj'). 



(E8) 



Finally, it should be noticed that p(co)q(co' ) is not a n ob servable in general, but the symmetrized product is Hermitian, 
and therefore, an observable. From Eqs. (E7) and (E8), we find 



(p(w)g(w') + «(w')p(w)> = 47rh 



n(w) 



Im{x(w)}o(w+w')- 



(E9) 



Equations (E7)-(E9) are general forms of the FDT. 



1. FDT for induced dipole fluctuations 

We apply the above general expressions of the FDT to dipole-dipole fluctuations, with the identifications 

p(t) Pi (t), 
q(t) Pj (t), 

X(t) OLij(t), 

<p{t) -+ Ej(t), 

where pi and pj are components of the dipole moment along directions i and j, respectively, Ej is the electric field 
along j at the position of the dipole, and ol^ is the (i, j) component of the polarizability tensor. The interaction 
Hamiltonian is H'(t) = —Ej(t)pj(t), where Ej is regarded as a time-dependent function and qj as an operator. The 
susceptibility acts in frequency space according to (Spi(u)) = aij(u)Ej(uj). With these substitutions, the FDT [Eqs. 
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(|E7k-(|E9|] takes the forms 



{pi{uj)pj{uj')) = Anh S(uj + u/) lm{aij(uj)} [rc(u;) + l], 
{pj{uj')pi{uj)) = Auk S(uj + uo r ) Im {aij(uj)} n(co>), 

^(PiMp^cjO+ft-^OPiM) = 47rfi5(o; + a; / )Ini{ai j (a;)} [n(w) + (ElO) 

From this formulation, the dissipation is found to be proportional to Im{a}. However, this leads to unphysical 
results in non-absorbing particles. For example, in spheres, non-vanishing values of Im{a} = (2k 3 /3)\a\ 2 arise from 



retardation corrections. Direct application of Eq. (Ell) to that case predicts an unphysical amount of absorption. In 
order to correct for this, we use a modified version of the FDT for spherical particles (a^- = SijQt), 

X - (piM^V) +Pj(u')Pi(u)) = 5(u + J) xH Sij [n(u) + i], (Ell) 
where x — Im{a} — (2& 3 /3)|a| 2 . The above derivation can be straightforwardly applied to magnetic dipole fluctuations. 

Appendix F: Fluctuation-dissipation theorem for electromagnetic fields 



Although it is possible to obtain an expression similar to Eq. (Ell) for the FDT corresponding to the fluctuations 
of the electromagnetic field following the above formalism [34], we present next an alternative derivation based upon 
the direct evaluation of electromagnetic quantum operator correlations. 

1. Electric field fluctuations 

The electric field quantum operator is defined as [42] 



E(r,o;) = E (+) (r,u;)+E ( } (r,u;), 



where 



,(+)/ \ • /27rfio; . ik . r 



E (r, lj) = i\l — —e Ka a Ka e zkr , (Fla) 



E ( - ) (r^) = -n/^e k ^ a e- k - (Fib) 



Here, k represents the wave vector, a is the polarization state, e^ ^/ is the polarization vector, and a k a (oJ ff ) is the 
photon annihilation (creation) operator, which acts over photon states as 

a k>' l nk ' a ^ = ( nk ' a + 1 )^k,k / ^,a' \(n + l)k,<x) , 

«kV |™k,<x) = n^aSkyS^a* \ (fl - l)k,a) • 

Using these expressions, we can calculate the correlation function of the electric field for normal ordering, 



E l ViX)®E l+; (r 2 ,a;")) = ]T (n k>CT |E l \r u u/) ® E l+j (r 2 ,u/') |n k>(7 ) (F2) 

k,cr,n 

= E(^3^3 / / dcouj2 K, CT |E ( ~ViV)®E (+) (r 2 VOK^>, 
where we assume a continuum of states and dQ denotes an element of solid angle in k space. The matrix element in 



Eq. (F2) reduces to 



(n k , a |E ( \r u u') (g)E (+) (r 2 ,cj // ) |n k , a ) = ^e k)(7 ®e k)ff e l ^' r n Ka 2n5(uj - uj')2ti5(uj - u"), 
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where r = ri — r 2 . Now, taking the average over all possible wave vectors and polarizations, we obtain 
^3^3 j dttJ2 J dujuj2 K^lE^^ri^O^E^^rs^'Olnk^) =fi J dtt^k 3 e Ka <g> e k>a e- ik " r n k><7 (5(a; / - u"). 

We use polarization vectors e p = and e s = (/>, which together with the wave vector k form an orthogonal basis set, 
such that 



e a ® e a = I 3 - k k. 



(F3) 



Moreover, we have 



and 



J dfte ±ikr = 2tt^ d(9sin((9)e ±i/crcosW = 2tt 



gifcr g— ikr 

ikr 



<ffi k (g) ke ±ik ' r 



2tt, 



V® V 



ikr 



Using all these results in Eq. ( F2 ) and assuming that the photons are at thermal equilibrium, we find 



E ( } (n, u)') <g> E (+) (r 2 , u)")\ = 47rfi(5(a/ - u") Im {0(n - r 2 , a;)} n(cj), 



(F4) 



where Q is the Green tensor defined in Eq. (C2) and n(uS) is the Bose-Einstein distribution function given by Eq. (E6). 



Following exactly the same procedure, we can calculate the correlation function of the electric field for antinormal 
ordering, 



E (+) (r 1;W ')®E ( ) (r 2;W ")}= £ <n k>ff |E ( \r 1 , J) ® E (+) (r 2 , a/') |n k)(T ) 

k,cr,n 

= E(^3 / / ^ {n Ka \t {+ \v l) uj')®t ( ~\v 2) uj'')\n^), 

which reduces to an expression similar to Eq. ( |F4[ ), except that n(uj) has to be substituted by n(uj) + 1: 
E (+) ( ri ,cj) ® E (_) (r 2 ,o/)\ = ^hd(uj - a/) lm{0(ri - r 2 ,u;)} [n(cj) + 1]. 



In a semiclassical framework such as the one used here (i.e., when the electromagnetic fields are classical magnitudes), 
the field fluctuations correspond to the symmetrized ordering, and therefore we have 

(E*(n,u;)<8>E(r2,a/)} = \ (^(n, ® E (+) (r 2 , a; 7 ) + E (+) (r 2 , cj) E (_) (n, a; 7 ) 
= 47tM(cj - a/) Im{£(ri - r 2 ,cj)} [n(u) + -]. 



2. Magnetic field fluctuations 



The Maxwell-Faraday law 
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c9t 



H(r,t) = V x E(r,t), 



allows us to obtain the magnetic field operator from the electric field [see Eqs. (Fl)]: 



H(r,cj) = H (+) (r,cj)+H ( } (r,u;), 
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where 



H (r, cj) = z\/ — — (k x e k?CT ) a k , a e 



(-)/_ A , / 27rfiw /i ^+ .-ik-r 



H (r, u) = — — (k x e k?CT ) a+ a e 



In way similar to the electrical field fluctuations, the correlation function of the magnetic field is given by 
(H ( ~VV)®H (+) (r2V')) = E K,.|H ( ~ViV)®H (+ V2V0K,.> 

k,cr,n 

= E(^^ I dn ^J dWOj2 M^~\r 1 ,w')®& + \r 2 ,u'')\n k , l ,) 

for normal ordering, and 

(H (+) (n , a/) ® H (_) (r 2 , w ") ) = E I H (+) (n , a/) ® H (_) (r 2 , a/') |n k><T > 

k,cr,n 

= E(^3 / d "E / dww2 <nk,a|H (+) (n, W , )®H ( " ) (r 2 ,w")|n k)ff > ) 

for antinormal ordering. From here, the calculation of the correlation functions is exactly the same as in the electric 
field. Incidentally, the sum over polarizations now becomes 



^(kx e Ka ^ (g) (k x e Ka ^j = ^ e k?a ® e kj<7 , 



so that we are again concerned with the same sum as in Eq. (F3), and therefore, the final result is the same as in the 
electric case: 

\ti,lj) <S> H^(r 2 , uj f )^ = 47rhS(uj — uo') lm {Q{r\ — r 2 ,co>)} n(co>), 
H (+) ( ri ,cj) H (_) (r2^')) = 4^hS(uj - a/) lm{0(n - r 2 ,u;)} [n(uj) + 1], 

^ (H (_) (n,cj) 0H (+) (r 2 ,^) +H (+) (r 2 ,cj) (g) H (_) (n, ^)) = 47rfiJ(o; - a/) Im {£(n - r 2 ,cj)} [n(u) + 

3. Electric- magnetic field fluctuations 

The correlation function of the electric-magnetic fields with normal ordering is 

jE ( " ) (r 1 , W ')®H (+) (r 2;W "))= E K-l E (_) (n, w') ® H (+) (r 2 , W ") |n k)<T > (F5) 

k,cr,n 

= E?i / d "E / dww2 K^Ie^^^O^h^^^'OK.). 

n ^ ' a 

From the Maxwell- Faraday law, we find 

H (+) (rV) = !vx E (+) (ry), (F6a) 
H 1 ' 1 ^') = -^V x E ( ' 1 (r > w ,, ) 1 (F6b) 
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Using Eq. (F6a), we can rewrite Eq. (F5) as 



E ( \ ri ,u;')®H (+ \r 2 ,u J ")) = ^V2xlj2j^^ J ' dQ Yl j ' ^ K )CT |E ( \t 1 ,cj') ® E (+) (r 2 , uj") |n k . ff ) } . 



and from here, we obtain 



E ( '(n.w)®^^,^)) = 4ttM(w - w') {V 2 x 0(n -r 2) w)} n(w). 



Likewise, for antinormal ordering, using Eq. (F6b), we find 



(F7) 



E (+) (r 1)W )(g>H ( \r 2 ,u')) = -47rM(w-w')4;Im{V2 x Q( ri -r 2) w)} [n(w) + 1]. (F8) 



4. Magnetic-electric field fluctuations 

The correlation function of the magnetic-electric fields for normal ordering can be written 
; H (_) (n , oj') ® E (+) (r 2 , a,")) = E |H (_) (ri,w') <g> E (+) (r 2 , a/') |n k>(T ) 

k,cr,n 



/ / dww2 K, CT |H ( ) (n,a; , )®E (+) (r 2 ,a;")|n k , <7 >. 



Then, using Eqs. (F6), we can transform this expression into 



1 „ v- V 



(27r) 3 c : 



JcM^ j dojtu 2 (n k , CT |E ( } (n, a/) ® E (+) (r 2 , w") |n k><T ) 



and from here, 



H ( \r 1 ,io) ® E (+) (r 2 ,w')) = -4ttM(w - w') 4 Im {Vi x £( r i " r 2> w )} 



(F9) 



Likewise, for antinormal ordering, we have 

r( + )/ A^-A(-) 

k,cr,n 



H W (r 1)W ')®E l J (r 2)W "))= £ (n k>ff | H (+) (r l5 a/) ® E ( ) (r 2 , to") \n K<r ) 

k,cr,n 

E(^3 / / <nk,a|H (+) (n,w , )®E ( " ) (r 2 ,w")|n k)ff > ) 



which, following the same procedure as in the previous sections and using Eqs. (F6), becomes 



H (+) (n,w)®E ( ) (r 2 ,o/)) =47rM(w-w') ^Im{Vi x0(n-r 2 ,w)} [n(w) + 1]. (F10) 



Finally, using Eqs. (F7)-(F10), we can write the symmetrized correlations of the electric- magnetic and the magnetic- 



electric fluctuations as 

X - (e { ~\t u uj) (g)H (+) (r 2 ,u/) +H (+) (r 2 ,o;) 0E H (r 1)W ')) = 47r&f(a; - a/) ^Im{Vi x 0(n -r 2 ,o;)} [n(cj) + 



\ (H ( } (n^) E (+) (r 2 ,^) + E (+) (r 2 ,cj) H ( ^n,"')) = -^HS(uj - a/) ^Im{Vi x Q( Tl - r 2 ,u;)} [n(uj) + J]. 
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Appendix G: Derivation of Eq. (4) 



We follow the notation and definitions introduced in Sec. II, and we refer to Fig. 1 there for the system under 
consideration. The net power absorbed by particle 1 can be written as [see Eqs. (1) and (2)] 



-( 



J — ( 



duoduo' 
duoduo' 



i{u-u')t i(jJ < (e+(J) \ia + {u') - (2fc 3 /3) Hw')| 2 ] -Ei(w)) 



(Gl) 



i(uj—uj')t 



2ur 

id (pf(uo f )g 12 (uo)p 2 {uo)) - (^(u/)piM) 



Now, using the expressions for the self-consistent fields and dipoles (see Sec. II) 

pi =^r x [rf + ai0i2pf] , 

where the deno mina tors D = I3 — ^12^2^21^1? -Di = I3 — ^1^12^2^21, and D2 = I3 — ^2^21^1^12 describe multiple 
scattering, Eq. (Gl) becomes 



Pi 



duoduo' 



(27T) 2 



i(o;— co )t 



Ef+ [ia+ - (2k 3 /3) |a| 2 ] D~ 1 Efj + (E* + [*a+ - (2fc 3 /3) |a| 2 ] D~ l g 12 a 2 E^ 

El+a+Qti - (2fc 3 /3) |a| 2 l -D _1 i?f 



E% + a+g+ 2 (D~ l ) + \ia+ - (2k 3 /Z) |a| 2 ] D^G^Eg 



We then carry out the matrix multiplications in this expression and apply the FDT. After some algebra, we find 



+ 2 



n 



2/ ' 



h + h 



J ° v=E,M 



Re< 










s 



Im { a v 2 C 
x ( n + 1 



g v a 2 A 



g v a 2 C 



g v a v 2 A^f + 0% A 



g v a 2 C 



v > n u v w v , 



2 , OL 2 A 







s 





Im {A} 



Re{C}+ -Rela^B 



lm{B} 



2 ' 



j 4 =- r^^fc 3 Yi xi 



2KI 2 



8 5 



+ 2 



a 2 



s 5 



n 
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7? f°° ( 



-A 2 |w„| 2 + g v ACu* v w v - g v ACu v w* v + A 2 \w„\ 2 



+ 7 



C 2 |w„| 2 - g v ACu* v w v + g v ACu v w* v - C 2 \w„\ 2 



2 2 

^3 43/ 



(afotf - o%af) - a% {a? af - af a™) (^C 3 - ^ 3 C) [> ( m + 

5 * rim {^ 



27' 



v=E,M 



\Au u \ 2 + g v AC*u v w v - g v AC*u v w* v - \Aw v \ 2 



\Cu v >\ 2 - g v A*Cu v ,w v > + g v A*Cu v >w* v , - \Cw u >\ 2 



\ {afaf - o$a?)* (a (C*) 3 - A* | A| 2 C*) 
v - 2 4 (a?a? - a™ a 2 M )* (a* |C| 2 C* - (A*) 3 6) } (ra 2 + , 



v=E,M 







2 




2 


1 


2" 




2 


5 


+ 2 


5 











a 



A^-C^ 



+ 



a, 



S S 



2 1 

+ 2 



ol\B 



n 2 + 



1 



v=E,M 

where the definition of the different variables is the same as in Sec. II. In this derivation, we have used the identity 

Im {xy 2 } = 2Re {xy} Im {y} + \y\ 2 Im {x} = 21m {xy} Re {y} — \y\ 2 Im {x} . 

Furthermore, the integral over uj' has been done using the Dirac 5 function of the FDT. After straightforward, lengthy 
algebraic manipulations, we can rearrange the above terms as 



V\ = #oi [no (oj) - ni (uj)} + H 02 [n (u) - n 2 (oj)} + H 2 i [n 2 (uj) - m (oj)} , 



(G2) 



where 



H 01 = - J2 [ duj4uxi 



+ 2 



+ 



-ImK} A^-g v C^ 



-Im 



V C^- 1" + Im {c^A 2 } ^ I" - Im {c^C 2 } I ^ " + 2^Re Im 



s I ' " 'Is " ^ ~ ' [ \s\ 2 



+Im {c^' A 2 } | ^ | 2 - Im {c^'C 2 } | y f - 2<fc,Re {qJ'Ac} Im j J + Re {a v 2 B} Im {5} 



02 
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- £ / 



v=E,M J ° 

duj4uj X Ulm{a^} 



A u v w v 

A g v C — 

s s 



+ 



a 2 



A^-g v C^ 
s s 



2 1 

+ - 



a 2 



,B 



A^-g„C^ 
s s 



v=E,M ' 

Finally, Eq. (4) can be readily obtained by reorganizing these expressions. 



+ Im \*i } \A^ - g v C^- 2 + ilm {a» 2 } 
I ) I s s 2 
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FIG. 8: (a,b) Heat transfer coefficient (HTC) as a function of particle distance d and temperature T for a gold-SiC dimer 
(solid curves) a nd a SiC-god dimer (broken curves). The HTC from the right particle (2) to the left particle (1) is defined as 
Vi/ST [see Eq. (G2)] with T\ = T = T and T 2 = T + 5T [see insets in (c) for sketches of the dimers]. (c,d) Power emanating 
from the right particle in the two dimer orientations under consideration. 



Appendix H: Symmetric and asymmetric heat transfers in inhomogeneous dimers 



We show in Fig. (8^a,b) the different heat transfer coefficients (HTCs, see Sec. Ill) exhibited by an inhomogeneous 
gold-SiC dimer when the hotter particle is either SiC (solid curves) or gold (dashed curves) and the transfer is for 
heat absorbed by the remaining gold or SiC particle, respectively. The power transfer to the cooler particle is nearly 
independent on whether the hotter particle is SiC or Au. This symmetry upon permutation of particle indices (1 <H> 2) 
is complete if we neglect multiple scattering between the particles [i.e., by setting u v — 1 and w v = in Eq. (4), since 
only and I^m terms contribute to the HTC], and indeed multiple scattering can be neglected in the clusters under 
consideration, as shown in Figs. 2-4. 

In contrast, the total heat lost by the hotter particle is strongly dependent on whether this is gold or SiC [see Fig. 
[8jc,d)]. Part of this heat is absorbed by the cooler particle, but the rest is radiated into the surrounding vacuum. At 
high temperatures above ~ 100 K, a hotter SiC particle produces larger radiation rates [Fig. [8jd)], rather independent 
of particle distance [Fig. |8jc)]. However, hotter gold is more capable of radiating at smaller temperatures [cf. vertical 
scales in Fig. |8jb,d)]. 



Appendix I: Total rate of heat loss and fraction of power exchange 



Heat exchange with the surrounding vacuum plays a leading role in the energy balance during radiative thermal- 
ization of a particle dimer, as shown in Fig. [9] We plot in Fig. [9ja) the power lost by the right particle when it 
is placed at a slightly higher temperature with respect to both the vacuum and the left particle (i.e., we represent 
—V2/ST for T = T\ = T and T 2 = T + ST). At high temperatures, this power is rather insensitive to the presence 
of the neighboring particle and no significant dependence on particle separation d is observed. This happens above 
T ~ 50 K in the inhomogeneous gold-SiC dimer and above T ~ 300 K in the homogeneous gold-gold dimer. However, 
the accompanying left particle has a strong influence on the power loss of the right particle at lower temperatures, 
specially in the homogeneous dimer. The power loss increases when the particles are placed closer together. We show 
next that this is partially explained by the effect of additional absorption by the left particle. 

Figure (9^b) shows the fraction of the power lost by the right particle that is absorbed by the left particle. This 
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FIG. 9: (a) Heat power lost by the right particle as a function of temperature for different particle separations. The right 
particle is slightly hotter than both the vacuum and the left particle, which are prepared at the same temperature, and the 
rate is normalized to the small temperature difference of the right particle (i.e., we represent —V2/8T for To = T\ — T and 
T2 = T + ST), (b) Fraction of the power lost by the right particle that ends up being absorbed by the left particle. 



fraction drops to small values at large temperatures, but it eventually approaches 100% at lower temperatures. This 
behavior is consistent with the distance dependence of the power lost by the right particle. Nearly full radiative 
heat transfer between two neighboring particles with negligible radiation into the surrounded vacuum is thus possible 
at sufficiently low temperatures. Specifically, in the gold-gold dimer this regime is already achieved at ~ 100 K for 
particles of radius R = 500 nm and a surface-to-surface separation of one radius (i.e., d = 1.5 fim). We observe that 
the temperature below which nearly 100% transfer between the dimer particles takes place decreases with increasing 
separation. 
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